home *** CD-ROM | disk | FTP | other *** search
/ Cream of the Crop 26 / Cream of the Crop 26.iso / os2 / octa209s.zip / octave-2.09 / libcruft / lapack / zlarfb.f < prev    next >
Text File  |  1996-07-19  |  19KB  |  610 lines
  1.       SUBROUTINE ZLARFB( SIDE, TRANS, DIRECT, STOREV, M, N, K, V, LDV,
  2.      $                   T, LDT, C, LDC, WORK, LDWORK )
  3. *
  4. *  -- LAPACK auxiliary routine (version 2.0) --
  5. *     Univ. of Tennessee, Univ. of California Berkeley, NAG Ltd.,
  6. *     Courant Institute, Argonne National Lab, and Rice University
  7. *     September 30, 1994
  8. *
  9. *     .. Scalar Arguments ..
  10.       CHARACTER          DIRECT, SIDE, STOREV, TRANS
  11.       INTEGER            K, LDC, LDT, LDV, LDWORK, M, N
  12. *     ..
  13. *     .. Array Arguments ..
  14.       COMPLEX*16         C( LDC, * ), T( LDT, * ), V( LDV, * ),
  15.      $                   WORK( LDWORK, * )
  16. *     ..
  17. *
  18. *  Purpose
  19. *  =======
  20. *
  21. *  ZLARFB applies a complex block reflector H or its transpose H' to a
  22. *  complex M-by-N matrix C, from either the left or the right.
  23. *
  24. *  Arguments
  25. *  =========
  26. *
  27. *  SIDE    (input) CHARACTER*1
  28. *          = 'L': apply H or H' from the Left
  29. *          = 'R': apply H or H' from the Right
  30. *
  31. *  TRANS   (input) CHARACTER*1
  32. *          = 'N': apply H (No transpose)
  33. *          = 'C': apply H' (Conjugate transpose)
  34. *
  35. *  DIRECT  (input) CHARACTER*1
  36. *          Indicates how H is formed from a product of elementary
  37. *          reflectors
  38. *          = 'F': H = H(1) H(2) . . . H(k) (Forward)
  39. *          = 'B': H = H(k) . . . H(2) H(1) (Backward)
  40. *
  41. *  STOREV  (input) CHARACTER*1
  42. *          Indicates how the vectors which define the elementary
  43. *          reflectors are stored:
  44. *          = 'C': Columnwise
  45. *          = 'R': Rowwise
  46. *
  47. *  M       (input) INTEGER
  48. *          The number of rows of the matrix C.
  49. *
  50. *  N       (input) INTEGER
  51. *          The number of columns of the matrix C.
  52. *
  53. *  K       (input) INTEGER
  54. *          The order of the matrix T (= the number of elementary
  55. *          reflectors whose product defines the block reflector).
  56. *
  57. *  V       (input) COMPLEX*16 array, dimension
  58. *                                (LDV,K) if STOREV = 'C'
  59. *                                (LDV,M) if STOREV = 'R' and SIDE = 'L'
  60. *                                (LDV,N) if STOREV = 'R' and SIDE = 'R'
  61. *          The matrix V. See further details.
  62. *
  63. *  LDV     (input) INTEGER
  64. *          The leading dimension of the array V.
  65. *          If STOREV = 'C' and SIDE = 'L', LDV >= max(1,M);
  66. *          if STOREV = 'C' and SIDE = 'R', LDV >= max(1,N);
  67. *          if STOREV = 'R', LDV >= K.
  68. *
  69. *  T       (input) COMPLEX*16 array, dimension (LDT,K)
  70. *          The triangular K-by-K matrix T in the representation of the
  71. *          block reflector.
  72. *
  73. *  LDT     (input) INTEGER
  74. *          The leading dimension of the array T. LDT >= K.
  75. *
  76. *  C       (input/output) COMPLEX*16 array, dimension (LDC,N)
  77. *          On entry, the M-by-N matrix C.
  78. *          On exit, C is overwritten by H*C or H'*C or C*H or C*H'.
  79. *
  80. *  LDC     (input) INTEGER
  81. *          The leading dimension of the array C. LDC >= max(1,M).
  82. *
  83. *  WORK    (workspace) COMPLEX*16 array, dimension (LDWORK,K)
  84. *
  85. *  LDWORK  (input) INTEGER
  86. *          The leading dimension of the array WORK.
  87. *          If SIDE = 'L', LDWORK >= max(1,N);
  88. *          if SIDE = 'R', LDWORK >= max(1,M).
  89. *
  90. *  =====================================================================
  91. *
  92. *     .. Parameters ..
  93.       COMPLEX*16         ONE
  94.       PARAMETER          ( ONE = ( 1.0D+0, 0.0D+0 ) )
  95. *     ..
  96. *     .. Local Scalars ..
  97.       CHARACTER          TRANST
  98.       INTEGER            I, J
  99. *     ..
  100. *     .. External Functions ..
  101.       LOGICAL            LSAME
  102.       EXTERNAL           LSAME
  103. *     ..
  104. *     .. External Subroutines ..
  105.       EXTERNAL           ZCOPY, ZGEMM, ZLACGV, ZTRMM
  106. *     ..
  107. *     .. Intrinsic Functions ..
  108.       INTRINSIC          DCONJG
  109. *     ..
  110. *     .. Executable Statements ..
  111. *
  112. *     Quick return if possible
  113. *
  114.       IF( M.LE.0 .OR. N.LE.0 )
  115.      $   RETURN
  116. *
  117.       IF( LSAME( TRANS, 'N' ) ) THEN
  118.          TRANST = 'C'
  119.       ELSE
  120.          TRANST = 'N'
  121.       END IF
  122. *
  123.       IF( LSAME( STOREV, 'C' ) ) THEN
  124. *
  125.          IF( LSAME( DIRECT, 'F' ) ) THEN
  126. *
  127. *           Let  V =  ( V1 )    (first K rows)
  128. *                     ( V2 )
  129. *           where  V1  is unit lower triangular.
  130. *
  131.             IF( LSAME( SIDE, 'L' ) ) THEN
  132. *
  133. *              Form  H * C  or  H' * C  where  C = ( C1 )
  134. *                                                  ( C2 )
  135. *
  136. *              W := C' * V  =  (C1'*V1 + C2'*V2)  (stored in WORK)
  137. *
  138. *              W := C1'
  139. *
  140.                DO 10 J = 1, K
  141.                   CALL ZCOPY( N, C( J, 1 ), LDC, WORK( 1, J ), 1 )
  142.                   CALL ZLACGV( N, WORK( 1, J ), 1 )
  143.    10          CONTINUE
  144. *
  145. *              W := W * V1
  146. *
  147.                CALL ZTRMM( 'Right', 'Lower', 'No transpose', 'Unit', N,
  148.      $                     K, ONE, V, LDV, WORK, LDWORK )
  149.                IF( M.GT.K ) THEN
  150. *
  151. *                 W := W + C2'*V2
  152. *
  153.                   CALL ZGEMM( 'Conjugate transpose', 'No transpose', N,
  154.      $                        K, M-K, ONE, C( K+1, 1 ), LDC,
  155.      $                        V( K+1, 1 ), LDV, ONE, WORK, LDWORK )
  156.                END IF
  157. *
  158. *              W := W * T'  or  W * T
  159. *
  160.                CALL ZTRMM( 'Right', 'Upper', TRANST, 'Non-unit', N, K,
  161.      $                     ONE, T, LDT, WORK, LDWORK )
  162. *
  163. *              C := C - V * W'
  164. *
  165.                IF( M.GT.K ) THEN
  166. *
  167. *                 C2 := C2 - V2 * W'
  168. *
  169.                   CALL ZGEMM( 'No transpose', 'Conjugate transpose',
  170.      $                        M-K, N, K, -ONE, V( K+1, 1 ), LDV, WORK,
  171.      $                        LDWORK, ONE, C( K+1, 1 ), LDC )
  172.                END IF
  173. *
  174. *              W := W * V1'
  175. *
  176.                CALL ZTRMM( 'Right', 'Lower', 'Conjugate transpose',
  177.      $                     'Unit', N, K, ONE, V, LDV, WORK, LDWORK )
  178. *
  179. *              C1 := C1 - W'
  180. *
  181.                DO 30 J = 1, K
  182.                   DO 20 I = 1, N
  183.                      C( J, I ) = C( J, I ) - DCONJG( WORK( I, J ) )
  184.    20             CONTINUE
  185.    30          CONTINUE
  186. *
  187.             ELSE IF( LSAME( SIDE, 'R' ) ) THEN
  188. *
  189. *              Form  C * H  or  C * H'  where  C = ( C1  C2 )
  190. *
  191. *              W := C * V  =  (C1*V1 + C2*V2)  (stored in WORK)
  192. *
  193. *              W := C1
  194. *
  195.                DO 40 J = 1, K
  196.                   CALL ZCOPY( M, C( 1, J ), 1, WORK( 1, J ), 1 )
  197.    40          CONTINUE
  198. *
  199. *              W := W * V1
  200. *
  201.                CALL ZTRMM( 'Right', 'Lower', 'No transpose', 'Unit', M,
  202.      $                     K, ONE, V, LDV, WORK, LDWORK )
  203.                IF( N.GT.K ) THEN
  204. *
  205. *                 W := W + C2 * V2
  206. *
  207.                   CALL ZGEMM( 'No transpose', 'No transpose', M, K, N-K,
  208.      $                        ONE, C( 1, K+1 ), LDC, V( K+1, 1 ), LDV,
  209.      $                        ONE, WORK, LDWORK )
  210.                END IF
  211. *
  212. *              W := W * T  or  W * T'
  213. *
  214.                CALL ZTRMM( 'Right', 'Upper', TRANS, 'Non-unit', M, K,
  215.      $                     ONE, T, LDT, WORK, LDWORK )
  216. *
  217. *              C := C - W * V'
  218. *
  219.                IF( N.GT.K ) THEN
  220. *
  221. *                 C2 := C2 - W * V2'
  222. *
  223.                   CALL ZGEMM( 'No transpose', 'Conjugate transpose', M,
  224.      $                        N-K, K, -ONE, WORK, LDWORK, V( K+1, 1 ),
  225.      $                        LDV, ONE, C( 1, K+1 ), LDC )
  226.                END IF
  227. *
  228. *              W := W * V1'
  229. *
  230.                CALL ZTRMM( 'Right', 'Lower', 'Conjugate transpose',
  231.      $                     'Unit', M, K, ONE, V, LDV, WORK, LDWORK )
  232. *
  233. *              C1 := C1 - W
  234. *
  235.                DO 60 J = 1, K
  236.                   DO 50 I = 1, M
  237.                      C( I, J ) = C( I, J ) - WORK( I, J )
  238.    50             CONTINUE
  239.    60          CONTINUE
  240.             END IF
  241. *
  242.          ELSE
  243. *
  244. *           Let  V =  ( V1 )
  245. *                     ( V2 )    (last K rows)
  246. *           where  V2  is unit upper triangular.
  247. *
  248.             IF( LSAME( SIDE, 'L' ) ) THEN
  249. *
  250. *              Form  H * C  or  H' * C  where  C = ( C1 )
  251. *                                                  ( C2 )
  252. *
  253. *              W := C' * V  =  (C1'*V1 + C2'*V2)  (stored in WORK)
  254. *
  255. *              W := C2'
  256. *
  257.                DO 70 J = 1, K
  258.                   CALL ZCOPY( N, C( M-K+J, 1 ), LDC, WORK( 1, J ), 1 )
  259.                   CALL ZLACGV( N, WORK( 1, J ), 1 )
  260.    70          CONTINUE
  261. *
  262. *              W := W * V2
  263. *
  264.                CALL ZTRMM( 'Right', 'Upper', 'No transpose', 'Unit', N,
  265.      $                     K, ONE, V( M-K+1, 1 ), LDV, WORK, LDWORK )
  266.                IF( M.GT.K ) THEN
  267. *
  268. *                 W := W + C1'*V1
  269. *
  270.                   CALL ZGEMM( 'Conjugate transpose', 'No transpose', N,
  271.      $                        K, M-K, ONE, C, LDC, V, LDV, ONE, WORK,
  272.      $                        LDWORK )
  273.                END IF
  274. *
  275. *              W := W * T'  or  W * T
  276. *
  277.                CALL ZTRMM( 'Right', 'Lower', TRANST, 'Non-unit', N, K,
  278.      $                     ONE, T, LDT, WORK, LDWORK )
  279. *
  280. *              C := C - V * W'
  281. *
  282.                IF( M.GT.K ) THEN
  283. *
  284. *                 C1 := C1 - V1 * W'
  285. *
  286.                   CALL ZGEMM( 'No transpose', 'Conjugate transpose',
  287.      $                        M-K, N, K, -ONE, V, LDV, WORK, LDWORK,
  288.      $                        ONE, C, LDC )
  289.                END IF
  290. *
  291. *              W := W * V2'
  292. *
  293.                CALL ZTRMM( 'Right', 'Upper', 'Conjugate transpose',
  294.      $                     'Unit', N, K, ONE, V( M-K+1, 1 ), LDV, WORK,
  295.      $                     LDWORK )
  296. *
  297. *              C2 := C2 - W'
  298. *
  299.                DO 90 J = 1, K
  300.                   DO 80 I = 1, N
  301.                      C( M-K+J, I ) = C( M-K+J, I ) -
  302.      $                               DCONJG( WORK( I, J ) )
  303.    80             CONTINUE
  304.    90          CONTINUE
  305. *
  306.             ELSE IF( LSAME( SIDE, 'R' ) ) THEN
  307. *
  308. *              Form  C * H  or  C * H'  where  C = ( C1  C2 )
  309. *
  310. *              W := C * V  =  (C1*V1 + C2*V2)  (stored in WORK)
  311. *
  312. *              W := C2
  313. *
  314.                DO 100 J = 1, K
  315.                   CALL ZCOPY( M, C( 1, N-K+J ), 1, WORK( 1, J ), 1 )
  316.   100          CONTINUE
  317. *
  318. *              W := W * V2
  319. *
  320.                CALL ZTRMM( 'Right', 'Upper', 'No transpose', 'Unit', M,
  321.      $                     K, ONE, V( N-K+1, 1 ), LDV, WORK, LDWORK )
  322.                IF( N.GT.K ) THEN
  323. *
  324. *                 W := W + C1 * V1
  325. *
  326.                   CALL ZGEMM( 'No transpose', 'No transpose', M, K, N-K,
  327.      $                        ONE, C, LDC, V, LDV, ONE, WORK, LDWORK )
  328.                END IF
  329. *
  330. *              W := W * T  or  W * T'
  331. *
  332.                CALL ZTRMM( 'Right', 'Lower', TRANS, 'Non-unit', M, K,
  333.      $                     ONE, T, LDT, WORK, LDWORK )
  334. *
  335. *              C := C - W * V'
  336. *
  337.                IF( N.GT.K ) THEN
  338. *
  339. *                 C1 := C1 - W * V1'
  340. *
  341.                   CALL ZGEMM( 'No transpose', 'Conjugate transpose', M,
  342.      $                        N-K, K, -ONE, WORK, LDWORK, V, LDV, ONE,
  343.      $                        C, LDC )
  344.                END IF
  345. *
  346. *              W := W * V2'
  347. *
  348.                CALL ZTRMM( 'Right', 'Upper', 'Conjugate transpose',
  349.      $                     'Unit', M, K, ONE, V( N-K+1, 1 ), LDV, WORK,
  350.      $                     LDWORK )
  351. *
  352. *              C2 := C2 - W
  353. *
  354.                DO 120 J = 1, K
  355.                   DO 110 I = 1, M
  356.                      C( I, N-K+J ) = C( I, N-K+J ) - WORK( I, J )
  357.   110             CONTINUE
  358.   120          CONTINUE
  359.             END IF
  360.          END IF
  361. *
  362.       ELSE IF( LSAME( STOREV, 'R' ) ) THEN
  363. *
  364.          IF( LSAME( DIRECT, 'F' ) ) THEN
  365. *
  366. *           Let  V =  ( V1  V2 )    (V1: first K columns)
  367. *           where  V1  is unit upper triangular.
  368. *
  369.             IF( LSAME( SIDE, 'L' ) ) THEN
  370. *
  371. *              Form  H * C  or  H' * C  where  C = ( C1 )
  372. *                                                  ( C2 )
  373. *
  374. *              W := C' * V'  =  (C1'*V1' + C2'*V2') (stored in WORK)
  375. *
  376. *              W := C1'
  377. *
  378.                DO 130 J = 1, K
  379.                   CALL ZCOPY( N, C( J, 1 ), LDC, WORK( 1, J ), 1 )
  380.                   CALL ZLACGV( N, WORK( 1, J ), 1 )
  381.   130          CONTINUE
  382. *
  383. *              W := W * V1'
  384. *
  385.                CALL ZTRMM( 'Right', 'Upper', 'Conjugate transpose',
  386.      $                     'Unit', N, K, ONE, V, LDV, WORK, LDWORK )
  387.                IF( M.GT.K ) THEN
  388. *
  389. *                 W := W + C2'*V2'
  390. *
  391.                   CALL ZGEMM( 'Conjugate transpose',
  392.      $                        'Conjugate transpose', N, K, M-K, ONE,
  393.      $                        C( K+1, 1 ), LDC, V( 1, K+1 ), LDV, ONE,
  394.      $                        WORK, LDWORK )
  395.                END IF
  396. *
  397. *              W := W * T'  or  W * T
  398. *
  399.                CALL ZTRMM( 'Right', 'Upper', TRANST, 'Non-unit', N, K,
  400.      $                     ONE, T, LDT, WORK, LDWORK )
  401. *
  402. *              C := C - V' * W'
  403. *
  404.                IF( M.GT.K ) THEN
  405. *
  406. *                 C2 := C2 - V2' * W'
  407. *
  408.                   CALL ZGEMM( 'Conjugate transpose',
  409.      $                        'Conjugate transpose', M-K, N, K, -ONE,
  410.      $                        V( 1, K+1 ), LDV, WORK, LDWORK, ONE,
  411.      $                        C( K+1, 1 ), LDC )
  412.                END IF
  413. *
  414. *              W := W * V1
  415. *
  416.                CALL ZTRMM( 'Right', 'Upper', 'No transpose', 'Unit', N,
  417.      $                     K, ONE, V, LDV, WORK, LDWORK )
  418. *
  419. *              C1 := C1 - W'
  420. *
  421.                DO 150 J = 1, K
  422.                   DO 140 I = 1, N
  423.                      C( J, I ) = C( J, I ) - DCONJG( WORK( I, J ) )
  424.   140             CONTINUE
  425.   150          CONTINUE
  426. *
  427.             ELSE IF( LSAME( SIDE, 'R' ) ) THEN
  428. *
  429. *              Form  C * H  or  C * H'  where  C = ( C1  C2 )
  430. *
  431. *              W := C * V'  =  (C1*V1' + C2*V2')  (stored in WORK)
  432. *
  433. *              W := C1
  434. *
  435.                DO 160 J = 1, K
  436.                   CALL ZCOPY( M, C( 1, J ), 1, WORK( 1, J ), 1 )
  437.   160          CONTINUE
  438. *
  439. *              W := W * V1'
  440. *
  441.                CALL ZTRMM( 'Right', 'Upper', 'Conjugate transpose',
  442.      $                     'Unit', M, K, ONE, V, LDV, WORK, LDWORK )
  443.                IF( N.GT.K ) THEN
  444. *
  445. *                 W := W + C2 * V2'
  446. *
  447.                   CALL ZGEMM( 'No transpose', 'Conjugate transpose', M,
  448.      $                        K, N-K, ONE, C( 1, K+1 ), LDC,
  449.      $                        V( 1, K+1 ), LDV, ONE, WORK, LDWORK )
  450.                END IF
  451. *
  452. *              W := W * T  or  W * T'
  453. *
  454.                CALL ZTRMM( 'Right', 'Upper', TRANS, 'Non-unit', M, K,
  455.      $                     ONE, T, LDT, WORK, LDWORK )
  456. *
  457. *              C := C - W * V
  458. *
  459.                IF( N.GT.K ) THEN
  460. *
  461. *                 C2 := C2 - W * V2
  462. *
  463.                   CALL ZGEMM( 'No transpose', 'No transpose', M, N-K, K,
  464.      $                        -ONE, WORK, LDWORK, V( 1, K+1 ), LDV, ONE,
  465.      $                        C( 1, K+1 ), LDC )
  466.                END IF
  467. *
  468. *              W := W * V1
  469. *
  470.                CALL ZTRMM( 'Right', 'Upper', 'No transpose', 'Unit', M,
  471.      $                     K, ONE, V, LDV, WORK, LDWORK )
  472. *
  473. *              C1 := C1 - W
  474. *
  475.                DO 180 J = 1, K
  476.                   DO 170 I = 1, M
  477.                      C( I, J ) = C( I, J ) - WORK( I, J )
  478.   170             CONTINUE
  479.   180          CONTINUE
  480. *
  481.             END IF
  482. *
  483.          ELSE
  484. *
  485. *           Let  V =  ( V1  V2 )    (V2: last K columns)
  486. *           where  V2  is unit lower triangular.
  487. *
  488.             IF( LSAME( SIDE, 'L' ) ) THEN
  489. *
  490. *              Form  H * C  or  H' * C  where  C = ( C1 )
  491. *                                                  ( C2 )
  492. *
  493. *              W := C' * V'  =  (C1'*V1' + C2'*V2') (stored in WORK)
  494. *
  495. *              W := C2'
  496. *
  497.                DO 190 J = 1, K
  498.                   CALL ZCOPY( N, C( M-K+J, 1 ), LDC, WORK( 1, J ), 1 )
  499.                   CALL ZLACGV( N, WORK( 1, J ), 1 )
  500.   190          CONTINUE
  501. *
  502. *              W := W * V2'
  503. *
  504.                CALL ZTRMM( 'Right', 'Lower', 'Conjugate transpose',
  505.      $                     'Unit', N, K, ONE, V( 1, M-K+1 ), LDV, WORK,
  506.      $                     LDWORK )
  507.                IF( M.GT.K ) THEN
  508. *
  509. *                 W := W + C1'*V1'
  510. *
  511.                   CALL ZGEMM( 'Conjugate transpose',
  512.      $                        'Conjugate transpose', N, K, M-K, ONE, C,
  513.      $                        LDC, V, LDV, ONE, WORK, LDWORK )
  514.                END IF
  515. *
  516. *              W := W * T'  or  W * T
  517. *
  518.                CALL ZTRMM( 'Right', 'Lower', TRANST, 'Non-unit', N, K,
  519.      $                     ONE, T, LDT, WORK, LDWORK )
  520. *
  521. *              C := C - V' * W'
  522. *
  523.                IF( M.GT.K ) THEN
  524. *
  525. *                 C1 := C1 - V1' * W'
  526. *
  527.                   CALL ZGEMM( 'Conjugate transpose',
  528.      $                        'Conjugate transpose', M-K, N, K, -ONE, V,
  529.      $                        LDV, WORK, LDWORK, ONE, C, LDC )
  530.                END IF
  531. *
  532. *              W := W * V2
  533. *
  534.                CALL ZTRMM( 'Right', 'Lower', 'No transpose', 'Unit', N,
  535.      $                     K, ONE, V( 1, M-K+1 ), LDV, WORK, LDWORK )
  536. *
  537. *              C2 := C2 - W'
  538. *
  539.                DO 210 J = 1, K
  540.                   DO 200 I = 1, N
  541.                      C( M-K+J, I ) = C( M-K+J, I ) -
  542.      $                               DCONJG( WORK( I, J ) )
  543.   200             CONTINUE
  544.   210          CONTINUE
  545. *
  546.             ELSE IF( LSAME( SIDE, 'R' ) ) THEN
  547. *
  548. *              Form  C * H  or  C * H'  where  C = ( C1  C2 )
  549. *
  550. *              W := C * V'  =  (C1*V1' + C2*V2')  (stored in WORK)
  551. *
  552. *              W := C2
  553. *
  554.                DO 220 J = 1, K
  555.                   CALL ZCOPY( M, C( 1, N-K+J ), 1, WORK( 1, J ), 1 )
  556.   220          CONTINUE
  557. *
  558. *              W := W * V2'
  559. *
  560.                CALL ZTRMM( 'Right', 'Lower', 'Conjugate transpose',
  561.      $                     'Unit', M, K, ONE, V( 1, N-K+1 ), LDV, WORK,
  562.      $                     LDWORK )
  563.                IF( N.GT.K ) THEN
  564. *
  565. *                 W := W + C1 * V1'
  566. *
  567.                   CALL ZGEMM( 'No transpose', 'Conjugate transpose', M,
  568.      $                        K, N-K, ONE, C, LDC, V, LDV, ONE, WORK,
  569.      $                        LDWORK )
  570.                END IF
  571. *
  572. *              W := W * T  or  W * T'
  573. *
  574.                CALL ZTRMM( 'Right', 'Lower', TRANS, 'Non-unit', M, K,
  575.      $                     ONE, T, LDT, WORK, LDWORK )
  576. *
  577. *              C := C - W * V
  578. *
  579.                IF( N.GT.K ) THEN
  580. *
  581. *                 C1 := C1 - W * V1
  582. *
  583.                   CALL ZGEMM( 'No transpose', 'No transpose', M, N-K, K,
  584.      $                        -ONE, WORK, LDWORK, V, LDV, ONE, C, LDC )
  585.                END IF
  586. *
  587. *              W := W * V2
  588. *
  589.                CALL ZTRMM( 'Right', 'Lower', 'No transpose', 'Unit', M,
  590.      $                     K, ONE, V( 1, N-K+1 ), LDV, WORK, LDWORK )
  591. *
  592. *              C1 := C1 - W
  593. *
  594.                DO 240 J = 1, K
  595.                   DO 230 I = 1, M
  596.                      C( I, N-K+J ) = C( I, N-K+J ) - WORK( I, J )
  597.   230             CONTINUE
  598.   240          CONTINUE
  599. *
  600.             END IF
  601. *
  602.          END IF
  603.       END IF
  604. *
  605.       RETURN
  606. *
  607. *     End of ZLARFB
  608. *
  609.       END
  610.